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ABSTRACT 

The viscous evolution of a thin disc around a central object is considered. Such discs 
are described by self-similar solutions in which either all or none of the inflowing 
mass accretes. An approximate solution for the partial accretion regime is constructed 
by employing a prescription recently introduced for nonlinear heat conduction. The 
solution is compared with numerical simulations demonstrating that the approximate 
solution describes the intermediate asymptotic stage for partially accreting discs. 



1 INTRODUCTION 

Accretion of gaseous material onto a central obj ect is a process at the heart of many diverse astrophysical phenomena, from 
quasars to the formation of planetary systems JPrank et al. 2002l l . In many cases accretion is via a geometrically thin disc 
( Shakura fc SunvaevlllQT^ : iLvnden-Beh fc Pringlelll974l ). 

The angular velocity of matter in the disc, Q., is keplerian (Ok = \/ GM,/r^ where G is the gravitational constant, M* 
is the mass of the accreting object and r is the cylindrical radial distance from the center of the star) throughout the disc 
except in the narrow boundary layer where it matches to the angular velocity of the central object SI, at the inner radius of 
the disc, rin. 

Accreting matter carries with it angular momentum whi ch is acquired by the central object if it is rotating slower than 
the keplerian angular velocity at the inner radius of the disc (jPringle fc Reeslll972h . The rate of material angular momentum 
flux by accretion is imat = Mr^Q where M is the accretion rate. Viscous stresses also lead to a flux of angular momentum 
Lvis ~ 2-Rr''^vY,dQ,/ dr where u is the kinematic turbulent viscosity and E is the surface mass density. In the steady state the 
mass accretion rate, M, and the total angular momentum flux, L — Lmat + ivis, in the disc are integration constants. For a 
keplerian disc the angular momentum balance leads to 



j.E= — 
■in 



M 



L 



VGMTr, 

Usually one writes the total angular momentum flux in units of material stress for keplerian flow at the inner radius 



: l3MrlnK{r in) 



and the equation takes the more familiar form 



(1) 



(2) 



(3) 



Within the disc the viscous and the material torques can have different ratios though the total angular momentum flux 
is a constant determined by the fastness parameter to, = fi,/f2K('"in). If the central object is rotating slowly, the angular 
velocity reaches a maximum value before matching the angular velocity of the star. At this point the viscous torque, as it 
depends on the radial angular velocity gradient, vanishes and the material torque is equal to the total angular momentum 
flux. If the boundary layer is narrow this allows one to estimate /3 ~ 1 for slowly rotating objects. 

If the disc is truncated beyond the stellar radius e.g. by the strong magnetic dipole field of the star, it is possible without 
stellar break-up, that the angular velocity of the star is much greater than the keplerian angular veloc ity at the inner r adius 
of the disc (w* > 1). In such cases the mass flux can diminish due to the centrifugal barrier ( lUariono v fc Sun vaevI I QTsh and 



the vi scous torque is negative i.e. the angular momentum is transferred from the star to the disc (/3 < 0). ISunvaev fc Shakura 



19771 ) gave a steady state solutio n for such " dead discs " in w hich mass flux (and so the mate rial torque) is zero 



The steady state solutions of Shakura fc Sunyacv ( 19731 ) and Sunvaev fc Shakura ( 1977) have time-dependent counter- 
parts as given by Pringlel (1974). These are self-similar solutions and extend to the origin, r = 0. In the first solution the mass 
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accre tion rate declines with a power-law in time while the angular momentum of the disc is constant (L = 0) (|Cannizzo et al 
1990l) . In the second solution the mass of the disc is constant (M — 0) but angular momentum of the disc increases with the 



viscous torque at the inner boundary (|Pringlelll99ll ). These extreme boundary conditions and the extension of the disc to 
r = in the solution manifest the necessity that no length scale should exist in order that a self-similar solution can be found. 

In reality, as the disc is not infinitely thin, matter vertically away from the disc plane can accrete even when matter 
is propelled at the midplane. Thus there should be a transition stage between thes e full accretion and full propeller stages 
during which the central object will accrete a fraction of the inflowing material ( Menou et al. 19991 : Eksi fc Kutlu 2O10l : 



Romanova et al. 2004 ; Ustvugova et al.l 20061) depending on the fastness para meter. Furthermore, it is expected that an 



episodic accretion stage I Spruit fc Taam 19931 : D'Angelo fc Spruit 2011 . 2012 ) will be realized near the transition stage 
separating the full accretion and full propeller stages. Such a stage can not be described as a succession of steady states 
as M is an integration constant and hence whatever the accretion rate in the disc should be the rate of accretion onto the 
star. The self-similar solutions because they correspond to extreme boundary conditions, either L = or M = 0, also can 
not accommodate the continuum of angular momentum flux per mass accretion rate, L/ M , to be realized in the continuous 
transition between the full propeller and full accretor stages. 

In this work we construct an approximate self-similar time dependent solution via a prescription recently used by lEksi 

I 2009h for self-similar evolution of temperature distribution in a semi-infinite rod with partia l heat insulation at one boundary. 
In the next section disc st ructure equa tions and the self-similar solutions of Pringlel (1974) is reviewed. In §3 the "nonlinear 
superposition" method of lEksi (|2009l ) is applied to these solutions to construct a self-similar approximate solution with 
arbitrary values of /3. In §4 we demonstrate the accuracy of the approximate analytical solution by comparing it with numerical 
results. In the final section we discuss some implications for neutron star soft X-ray transients. 



2 EVOLUTION OF THE DISC 

We consider the evolution of a non-self-gravitating viscous thin disc in the gravitational field of a central mass M« . 



2.1 Thin Disc Equations 

The conservation of mass in cylindrical coordinates reads 

9E _ 1 dM 
dt 27rr dr 

where E is the surface mass density and 

M = 27rrE (-Vr) 

is the mass flow rate, Vr being the radial velocity. The conservation of angular momentum reads 



where Q. is the angular velocity of the fluid and 



or 



(4) 



(5) 



(6) 



(7) 



is the angular momentum flux (torque). The first term is the material torque and the second term is the viscous torque where 
u is the turbulent kinematic viscosity. 

The form of the conservation equations @ and ^ make it clear that in the steady state the accretion rate M and 
angular momentum flux L are two integrals of motion. Assuming the matter in the disc is in Keplerian orbits, Q. — Qk where 
Ok = \/ GM,/r'-^, Equation ((Gjl with the help of Equation Q reduces to 



or \ 

This then can be plugged into Equation (|4]) to give the diffusion equation 



as 

dt 



3d_ 

r dr 



1/2 d /■ i/2\ 

r Yr V^' ) 



(8) 



(9) 



for the evolution of the surface mass density. The equation is linear if u does not depend on E. In general the viscosity depends 
on E and the equation is nonlinear. 

The rest of the disc structure equations are as follows ([Frank et al ]|2002[ ) : We emplov lShakura fc Sunvae^ (| 19731 ) turbulent 
viscosity prescription v = aCsH where a ~ 0.01 — 0.1 is the dimensionless viscosity parameter, Cs = \/ P/p is the sound speed, 
H = Cs/^K is the disc scale height and P is the pressure. The density is then given by p = E/_ff. Assuming the heat 
produced by viscous processes are radiated locally, the energy balance can be written as 4crsBr*/3Efi: = fi^Efil where T 
is the temperature at the disc mid plane and k is the Rosseland mean opacity. These equations are then supported by the 
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equation of state P = P{p, T) and a prescription for the opacity k, = fi:(p, T). The pressure has contributions from both the gas 
pressure, Pgas = pk^T /fmip where p, is the mean molecular weight, and radiation pressure Prad = 4crsBr^/3c. We assume the 
gas pressure dominates the radiation pressure throughout the disc. Opacity is in general of the form k = kqp'^T'' where kq, a 
and b are constants depending on the dominating opacity regime (e.g. a = 6 = and kq ~ 0.4 cm^ g~^ for electron scattering). 
These algebraic disc structure equations can be solved among themselves first to obtain viscosity (and other variables) in 
terms of r and E in the form 

u = Cr^E" (10) 

where C, p and q are constants determined by the dominant opacity regime and pressure (see e.g. Cannizzo et al. ( 1990l ): 
1 20091 )). This then can be plugged in the diffusion equation (|9} to obtain an equation containin g E only. For discs 



Ertan et al 



in which gas pressure and electron scattering opacity dominates p — 1 and g = 2/3 ( Cannizzo et al. 1990l ). In comparing the 
analytical solutions with the numerical solutions we assume this regime to prevail throughout the disc though it is well known 
that different opacity and pressure regimes prevail at different locations of the disc depending on temperature T{r,t). 



2.2 Pringle Solutions 



In order to apply self-similarity methods we render the equation dimensionless via R — r/ro, r — t/to and a = E/Eo, and 
define uo ~ Ct-q'Sq. Choosing to = 4ro/3j^o we obtain 



da 
dr 



4 d 
RdR 



R 



1/2 



d_ 

OR 



R 



l/2^,+ l^ 



(11) 



Self-simi larity methods ( ZeTDojd^h&Raizei [l967|}_ pro vide two solution s of this equation as fir st found by IPringlel (|l974l ) 
(see also iLvnden-Bell fc Pringlel j 19741 ) : iTanakal {mm : iMineshigel l|l99ll ): iMineshige et al.1 (|l993l )). The first solution which 
we call the "full accretor solution" , for reasons to be explained in the following subsection, is 



a{R,T) 



where 



(1 + r) 5q-2p + 4 IfK^l _|_ T 



(1 + ^)^^'''^>^P(- 4(P-1S+.) ) 



1- kp,q[R{l+r) =<!~2P+ 



2 E_1 1/9 

<! + l 



if g / 0, 
if g = 



(4g-2p + 4)(5g-2p + 4)' 
The second solution which we call the "full propeller solution" is 



a{R,T) = 



(1 + r) 4<j-2p+4 ( 



4q-2p + 4 



P+1/2 
■J+1 



(1 + r)2^p-2) R-P' 



exp 



V 4(p-2)2(l + r); ' 



R{1+t) 4,-2p + 4 



5 P + 1/2 
2 q + 1 



1/9 



if g / 0, 
if g = 0. 



(12) 



(13) 



(14) 



Note that Equation ((9)1 is symmetric under translations in time and hence the solutions can also be shifted in time. We have 
exploited this to write the original solutions shifted as t — >■ r-|- 1. Note also that the solutions for the linear cases (g — 0) can 
be found from the nonlinear cas es by referring lim5_»o(l + Aq)^^'' = exp A. In the linear solutions the p — 2 case should be 
handled separately as is done by Pringlel ( 1974 ) . We do not write the solution for this case as it would not contribute to the 
main motivation of the paper. 

In the non-dimensionalization process we have employed four quantities, Eq, tq, to and vq, related by two equations 
= Ct-qEq and to = 4ro/3i^o- This means we are free to attribute any value for two of the variables e.g. Eq and tq. We fix 
these two quantities in terms of the initial mass and angular momentum of the disc (see Ertan et al. 20091 ). 



2.3 Properties of solutions 



The Pringle solutions given in Equations (|12|l and (|14p describe the evolution of thin discs with free outer boundaries and 
thus are more appropriate for discs that are not tru ncated e.g. by tidal torques. The physical meaning of the first solution 
(full accretor) was clarified bv lCannizzo et aL ( 199C ) who modele d the evolution of a disc formed by tidally disrupted star 
around a supermassive black h ole. S oon afterwards I Pringlel (|l99ll ) used the second solution for des cribing the evolu t ion o f 



circumbinary discs. iMineshige (1991) used the first solution for the evolution of dwarf nova discs. Mineshige et al. ( 19931 ) 



used it for the evolutio n of the putative disc around the neutron star assumed to form in SN 87A. Since the original work of 
Chatteriee et al.1 (|2000l ) the first solution has been employed for describing fallback discs (e.g. iPerna et al 



20001: Ertan et al 



2OO9I ) around voung neutron stars. [Eksi & Alpar (2003) employed both two solutions, the first one for the accretion and the 
second one for the propeller stages of the system, respectively. 
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Figure 1. Comparison of the analytical and numerical solutions of the diffusion equation (|9)l for the full accretor (panel a) and full 
propeller (panel b) inner boundary conditions, S = and 9(r^/^!^S)/9r = at r = rin, respectively. The filled circles show the numerical 
solution which starts from a Gaussian surface density distribution. On the left panel the solid line shows the approximate analytical 
solution given in Equation II30I I for n = 1. The dashed line shows the exact Pringle solution given in Equation I I19II . 



The properties of the solutions given in Equations (|12|l and (|14|l are understood and classified in terms of the time 
evolution they result for the mass 



Md = / E • 2'Krdr (15) 
Jo 

and the angular momentum 

Ld = / r^^xS • 2-Krdr, (16) 



of the disc. Here rout is the freely expanding outer boundary of the disc defined as the location where the square bracket 
terms in the solutions given in Equations (|12p and (|14p vanish. The lower limit of the integrals is zero as the solutions extend 
to the origin. A real disc will have a finite inner radius and so the integrals then can provide the correct parameters as long 
as this inner radius is much smaller than the outer radius. 



2.3.1 Full Accretor Solution 



According to the first solution given in Equation (|12p . the outer boundary of the disc evolves as rout oc (1 + tY/'-''^ 2p+4)^ 
The mass of the disc decays with Md oc (1 + -7-)^i/(59-2p+4) that mass fiux out of the disc is 



Md 



Mo(l + r) ^ 5<j-2p+4^ 



(17) 



where Mo = Svri'oEo is the initial mass flow rate in the disc. In this solution the angular momentum of the disc is constant 
(Ld = 0): the outer radius of the disc expands at the exact rate to take up the angular momentum lost from the inner disc 
by t he accreting ma ss. 

Cannizzo et al. 1 199Ci ) showed that this solution describes the evolution of the disc with the absorbing inner boundary 
condition E(rin,t) = at the inner radius of the disc, Hn. In this solution cr oc for most part of the disc except 

near the outer boundary which the te rms in square brackets ha ndle. As shown in the next section this also is the radial 
profile in the steady-state counterpart ( Shakura fc Sunvaev 19731 ) away from the inner boundary of the disc. The sel f-similar 
solut i ons are thus said to describe a quasi-static regime i.e. successive steady-states with decreasing accretion rates (jPringle 
1974 : IChan et al.ll2005l ). According to this solution the disc extends to the origin and a blows out t here. In the numerical 
solution (J overturns trying to match the inner boundary condition. As first remarked by Ertan et al. ( 20091 ) the time shifted 
solutions where ro and Eo are fixed by referring to the initial mass and angular momentum in the disc gives a better match 
to the numerical solution. The numerical solution settles to the self-similar solut ion in a few diffusive timescales. As such 
self-similar solutions describe the intermediate asymptotic regime (|Barenblattlll996l ') as they prevail after the initial conditions 
are forgotten and before the disc vanishes completely. 



2.3.2 Full Propeller Solution 

According to second solution given in Equation (|14p but shifted in time the outer radius evolves as rout oc (1 + T-)2/(4q-2p+4) 
The mass of the disc is constant (Md = 0) and the angular momentum of the disc increases as Ld oc (1 -f- -7-)i/(*'J-2p+4) w^qx. 
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the angular momentum flux evolves as 



id = 1/0(1 + r)4<j-2p+4-\ (18) 

where io = Mo^/GMro is the initial angular momentum flux. Pringlel (1991) showed that this describes the intermediate 
asymptotic regime of the evolution of a disc with the bound ary cond i tion d {vT,r^^^) / dr = at r = nn which allows no mass 



flux out of the d isc from the inner boundar y (see Eqn.(|5))). iPringld ([198l|) mentioned that this solution could describe the 



propeller regime l lUarionov fc Sunvaev|[l975 ) of a disc in which the centrifugal barrier set by a rapidly rotating magnetosphere 
prohibits accretion. In this solution radial profile is a oc 

^-(p+l/2)/{,+l) ^1^.^!^ 

again, e xtends to the origin and blows up. 
This also is the radial profile in the steady-state counterpart (jSunvaev fc Shakuralll977r ) as shown in the next section. The 
torque given in Equation ((TJ consists of two parts, a material and a viscous torque. As there is no matter flux out of the 
disc, in this solution, the angular momentum flux above is due to the viscous torque only. In order to take up the angular 
momentum transfered to the disc by this torque, the outer radius of the disc expands more rapidly in this solution than it 
does in the full accretor solution (see Figure [T]). 



3 AN APPROXIMATE SOLUTION WITH GENERAL BOUNDARY CONDITIONS 

The Pringle solutions correspond to two extreme stages: either all matter reaching the inner radius of the disc accretes while 
the angular momentum of the disc is constant (flrst solution) or none accretes (mass of the disc is constant) while the angular 
momentum of the disc increases with the viscous torque from the inner boundary (second solution). 

It would be desirable to have a solution which could accommodate the continuous range of L/Af. In the linear case [q = 0) 
such a solution can be constructed by superposing the solutions; namely one can multiply the solutions with appropriate 
constants and adds them up. In the nonlinear case superposition can not be applied. 

The heat equation dtu = d^X^xU with power-law heat conduction coefflcient (x = Xou'^) is mathematically a special 
case (obtained for 2p = 3q) of the diffusion equation describing the evol ution of the disc. There a re self-similar solutions 



describing the evolution of temperature on a semi-inflnite rod x = [0, oo) (jZel'Dovich fc Raizeilll967l ) which are then special 
cases of the Pringle solutions. The flrst solution is a special case of the full accretion solution and accommodates u = 
at X = and so all heat inflowing there is absorbed. The other solution is a special case of the full propeller solution and 



acco mmodates du/dx = at a; = which allows no heat flux through the inner boundary (perfect insulator). Recently, ( Ekgil 



[20091 ) constructed an approximate self-similar solution for the partial insulator boundary condition at a; = 0. In this section we 
are going to employ the same prescription to the Pringle solutions to obtain an approximate solution for the partial accretion 
case. 



3.1 Solutions in terms of fluxes 



5<j-2p+4 a_g 



As a flrst step, we write the flrst and second solutions in terms of mass and angular momentum flux, respectively. It is possible 
to write the flrst solution given in Equation (|12l) in terms of the dimensionless mass flux m = M/Mq = (1 + t)~^~ 

[ (ii-f rh) exp {- 



^p,q l+r 



R~ 



4(p-2)^(l + r) 



if g / 0, 
if g = 



(19) 



The second solution given in Equation ([14]), in terms of dimensionless angular momentum flux I = L/Lq — {1 + t) 49-2P+4 
can be written as 



\R-P{-iR-^^'')] exp 



^-11 



4(p-2)2(1+t) 



if g / 0, 
if g = 



(20) 



Note that the viscous torque is negative when this solution is valid so that —£ is actually positive. 

We have eliminated most occurrences of r in favor of mass flux in the full accretor solution and angular momentum flux 
in the full propeller solution. Note the "dual form" of the solutions (rh — >■ —£R~^^^) which would not exist if we insisted in 
eliminating all occurrences of r in favor of m and £, respectively. This leads to another useful property: The solutions as they 
are given in Equations H9\) and i2(J\) satisfy the diffusion equation even when m and £ are defined as 

1 



and 



= mo(l + tY 



^0(1 + r) 



M = 1 + 



A : 



5g - 2p + 4 
1 



(21) 



(22) 



4g - 2p + 4 

i.e. when they are multiplied by some arbitrary constants rho and £0. This, in the linear cases (q — Q), is equivalent to 
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multiplying the solution itself with a constant as the flux is then a multiplier of the solution. This property, indispensable for 
the superposition procedure, would not be valid if we insisted in eliminating all occurrences of r in favor of the fluxes. 



3.2 Comparing with the steady state solution 

The steady-state solution of E can be found from Equations ^ and (|1U|I . In dimensionless form this reads 

a= [i?-'' (m-ii?-^/')]^ . (23) 

Here m and £ are integration constants. The solutions in terms of fluxes, given in equations (|19|) and (|20|) . are equivalent to 
the original solutions given in Equations (|I2|) and (|I4|) . respectively. The meaning of the self-similar solutions when written 
in the latter form is much transparent showing clearly that the self-similar solutions are special solutions in the sense that 
they each have a single integration constant (in space) whereas any ratio between m and £ is possible in the the steady state 
solution given in Equation (|23p . The steady state counterparts of the solutions given in Equations l|19|l and (|20p would be 
a = [R-P{m)]^^^''+'^^ and a = [_R"''(-i7?"^/^)]^''(''+^' , respectively. The first of these is the full accretor solution with no 
viscous torque and inner radius at the origin. The other limiting case, a = \R -''(- £R-^^^)]^^^''+'^\ obtained for m = from 
Equation (|23l) . is known as the "dead disc" solution of ISunvaev fc Shakural (|l977l ). Such a disc is expected to form in the 
propeller regime (jlllarionov fc Sunvaev 1 19751 ) though it would not be steady. 



3.3 Nonlinear "superposition" 

The idea of going back from the limiting cases, to the general solution given in Equation (f23)) suggests a "superposition" 
principle 

ll/(<3-H) 

a^[R--(rh-£R-^^')y'''^'' (24) 



= R-P m 



\r-p 



-£R-^/^ 



1/(9 + 1) 



which reduces to ordinary summation in the linearity limit (g = 0). Using this procedure, the solutions given in Equations 
H19p and (|20p can be "superposed" as 



iR-Pim-lR-^^^)] exp 



P'l l + T 



R-P{m- £R'^''') 



"4(p-2)'-i(l+r); ' 



if g / 0, 
if g = 0. 



(25) 



Clearly, = case will give the "full accretor" and mo = will give the "full propeller" solution given in Equations (|19p and 
()20p . respectively. 

If we plug in the expression given in Equation H25|) into the diffusion equation pip we find that the condition for satisfying 
it reduces to 



q+1 



{l + r)-\f{R,T) 



+ ^1 ] R m 



m(l + r) £{1 + t) 



-2]RM+{'^^U^ + \]P 



0. 



(26) 



m / \ m £ J \ £ 

where f{R,T) is some complicated function not very important for our discussion. The first and the third terms in square 
brackets vanish upon using Equations (|2ip and (|22p . respectively. The second term does not vanish and the condition for 
Equation (|25p to satisfy the diffusion equation pip reduces to 



J^-!^RfiR,r) 
q+ll+T 



- 



(27) 



where we have used 2 — — A ; 



Thus, the expression given in Equation (|25p is not an exact solution unless mo = or 



£o = (reducing to one of the already known solutions given in Equations p9p and (|2Up . respectively), or when kp^q — which 
is accomplished when g = i.e. the linearity limit. It thus appears that we have not gained much with the unified expression 
given in Equation (|25p as it is an exact solution only in the already known cases. Yet, we will show in the following that 
the unified expression given in Equation (|25p not only illuminates the meaning and complementarity of the known solutions 
but also motivates a very accurate approximate solution describing the partial accretion regime becoming exact for r — >■ oo, 
practically achieved in a few viscous time-scales. 



3.4 Why not an exact solution? 

Could there be an exact solution of the diffusion equation which could accommodate boundary conditions in which angular 
momentum flux per u nit mass accretion rate could take any value? In order that a self-similar solution is found no scale should 
exist in the problem ( Barenblatt 199^). The expression given in Equation (|25p or any expression with both mass accretion 
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rate and angular momentum flux having finite values can not be an exact solution because of the presence of the length scale 
I = E(0)/S'(0) oc {L/Mf/GM creeping into the problem through the boundary condition and destroying the very condition 
for the existence of a self-similar solution. Accordingly, the exact solutions given in Equations (|19[) and \2Q\ are found either 
when i = (L = 0) or when / — >■ oo (M = 0). Any attempt to find a solution with a finite inn er radius also would suffer 
from the same problem as this also would introduce a length scale into the problem l|Spruit fc Ta am 2001) except for the 
case 2p = 3g where the equation is symmetric under translations so that the origin can be shifted. We stress that a solution 
obtained by the non-linear superposition of self-similar solutions is not equivalent to shifting the origin of the coordinates. 

The superposition procedure in the linear case also introduces a length scale into the problem. Yet the superposed 
expression = case in Equation (|25|) ') is an exact solution. The presence of a scale does not hinder the exact solution in 
the linear case because, as hinted in Equation H27|) . the length scale is coupled to the "non-linearity parameter" q: vanishing 
of either renders the expression given in Equation exact. 



3.5 A more useful "solution" 



Note that m and I have different time dependencies, as given by Equations H21|) and (|22|l , respectively. This implies the inner 
radius of the disc in the "solution" given in Equation (|25|) is time-dependent. Note, however, that the generalization to the 
expression given in Equation (|25|) is attained from solutions with zero inner radius. For example, the solution with no viscous 
torque but finite inner radius would be 



(28) 



where iiin = rin/ro (|Pringldll98ll ). rather than a = [R ^rh] 9+1 . In order to accommodate cases with finite angular momentum 



fiux we refer to Equation (|2]) defining the dimensionless torque /?. This in dimensionless form can be written as 
(recall Eqn. ([2])), and write the "solution" given in Equation (I25|l in the form 



\R-Pm{l- P{Ri^/R) 



1/2' 



l + T \ 



R-^rh{l - l3{Rin/R) 



l/2^^ 9+1 



1/9 



exp 



if g / 0, 
if g = 



(29) 



(30) 



"4(p-2)^(l + r) 

which accommodates the no- viscous-torque inner boundary condition for j3 — 1. The left panel of Figure [1] shows this 
approximate solution and the exact Pringle solution given in Equation (|19[1 together with the numerical solution obtained 
with the boundary condition z/E(rin,t) = 0. 

The numerical solution (see the next section for the details of the numerical method) is started from a Gaussian initial 
profile. In comparing the solutions we have calculated the mass and angular momentum in this Gaussian surface density 



distrib utio n and used thes e values for the analytical solution. No fitting procedure is applied. As first shown bv lCannizzo et al 



( 1990h and lPrindeJ |l99ll ) ,he figures demonstrate that once the initial conditions are forgotten in a few viscous timescales, 
the evolution of the disc proceeds self-similarly. It is instructive to compare the exact Pringle solution given in Equation (1191) 
and the approximate solution given in Equation (|30|) : The approximate solution describes the inner region of the disc much 
better though they are both satisfactory for the intermediate and outer regions. The ma ss flow rate i n the numerical and 
analytical solutions, however, are the same (see for example the right panel of Figure 1 in Ertan et al. 1 20091 ')'). In reality, in 
order that all inflowing mass accretes, the boundary condition E(rin,f) = must be satisfied. The exact self-similar solution 
overestimates the surface mass density near the inner region such that this extra matter at the inner region compensates 
the reduction in the mass fiow rate due to reduced density gradient. The right panel of Figure [l] stands for the full propeller 
solution. In this case we see that the exact solution given in Equation (|20p does very fine except that it extends to the origin 
though real discs have finite radius. The outer radius of the disc, in order to take up the angular momentum imparted from 
the inner boundary, moves out more rapidly in this case. 

The radial profile of the full propeller solution given in Equation (|20p . a ex. R~^^'^^^^''^^'''^^\ is steeper than the radial 
profile of the full accretor solution given in Equation ([T9|, a oc R-p/^i+^'> (jSpruit fc Taamlll993l ;). This leads to enhanced 
accretion in transitions from the propeller to the accreting regimes. This then will lead to the inwards movement of the 
inner disc radius rin l e ading to cyclic behavior in systems at the brink of transitions between accretor and propeller regimes 
I Spruit fc Taamlll993l : [D'Angelo fc Spruitll2012l V 



4 RESULTS: A SOLUTION FOR PARTIAL ACCRETION 

In this section we demonstrate that the approximate solution given in Equation pO[) describes discs with partial accretion 
accurately. 
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Figure 2. The left panel shows the surface mass density profile for dises with partial accretion. The cases for / = 0.01, 0.1, 0.5 and 1 
are shown. The right panel shows the relation between the fraction of mass flux, /, that traverses the inner boundary of the disc and the 
dimensionless torque at the inner boundary, /3 as given in Equation 1 140 I I. The data points are the f} values obtained by fitting the surface 
mass density given in Equation II30II to the numerical solution for a value of /. In the Figure we have not applied any fitting procedure, 
but used 5 = Ax/x-^ = 0.1576 from our simulation parameters. 



4.1 Numerical Method 

In solving the dif fusion equation (|9l) nu merically we first transform to tfie "specific angular momentum coordinate" x = r^^^ 
and use S = xTj ( Bath fc Pringlelll981f) . The equation becomes 

where i/ — Cx^''~'' S"^ . With this notation Equation ^ becomes 

M = 3n-^(iyS). (32) 
ox 

These equations are discretized and the new value of S at each grid point is found by 

Sr" = S, + + - 2(^.5),] . (33) 

(XiAx)^ 

The discretization of Eqn. (1321) ) gives accretion rate at each grid point as 

M = 3.i^^^H^^. (34) 
Ax 

We have used imax = 2000 grids Xi — xq + iAx where xo — r]^'^ . We have taken ro = 10* cm and ri^^^ — 10^^ cm. 
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4.2 Inner Boundary Condition for Partial Accretion 

The mass flux at the inner boundary of the disc (i — 0) is 

Aa; ^ ' 

Accordingly, the maximum mass flux is obtained for {1^3)0 ~ 0, the full accretor boundary condition, as 

Mmax = 3^^^. (36) 

We have used the boundary condition (i^5')o ~ for the full accretion regime (see left panel in Figure [T]). 

According to Equation (|35p again, A/* = if (1^5)0 = (j^5)i and we have employed this boundary condition for the full 
propeller regime (see right panel in Figure [T| . If a fraction / of Mmax accretes 

J-™ max = 37r (67) 

is satisfled which implies 

{uS)o = il' .f){iyS)i. (38) 

In the following, we employ this boundary condition for simulating the evolution of discs in which a fraction / of the inflowing 
mass flux traverses the inner boundary and the rest 1 — / remains in the disc. Note that f — 1 and f — cases correspond 
to the well know full accretor and full propeller boundary conditions, respectively. 
With X = r^l'^ and S = Ea; Equation ((3]) can be written as 

z/S= ^(a;-/3xi„) (39) 

Finding [yS)^ and {yS)\ from this equation and plugging into Equation (f38)l gives 

1-/3 

where 

^ _ As _ Ar 

which is the relative width of grids in the speciflc angular momentum space. For our simulation parameters 5 = 0.1576. 

We have performed simulations of disc evolution for a range of / values from 0.01 to 1 keeping the value of / constant 
throughout the simultation. For each value of / we have seen that the self-similar evolution of the disc is well described by 
the approximate solution given in Equation ((5n)) . for a certain value of /3 found by solving /3 from Equation (gD]): 

/? = 1 - i^<5. (42) 

The left panel of Figure [2] shows the surface density proflle in the inner region of the disc at t = 0.1 years for / = 0.01, 0.10, 
0.50 and 1.00. The corresponding /3 values are also written in the Figure. The right panel of Figure [2] shows the relation given 
in Equation H40p together with data points from our numerical simulations. The Figure also shows results for other values of 
(5 from 0.01 to 1. 

The relative width given in Equation (|40|) is determined by the number of grids we have chosen. It thus appears as a 
numerical quantity that would vanish if we could do numerical computing with inflnite number of grids. At the inner boundary 
of the disc the matter rotating with keplerian angular velocity has to come to corotation with the star over some region of 
narrow width, the boundary layer. In our numerical solution of the diffusion equation ((Oj in which the Keplerian rotation is 
built in, non-keplerian rotation in the boundary layer can not be modeled and it is implicitly assumed that this transition 
takes place at the innermost grid and is unresolved. For real discs which do not have grids, 5 has the physical meaning as 
the relative width of the region over which matter in the disc is brought into corotation with the star. Equation (|40|) then 
tells that the relative width of this transition region is determined by the fraction of matter that could traverse it and the 
dimensionless torque. The relative width is determined by the stresses that bring the matter in the disc to co-rotation with 
the star. It is not necessarily the viscous torques. In the case of magnetized stars truncating the disc beyond the stellar surface 
this could be the magnetic stress applied by the rotating magnetosphere of the star. 



/ (40) 

(41) 



4.3 Partial Accretion Stage Away From Equilibrium 

In the previous subsection the value of / was kept constant throughout the simulation. It was found that the accretion rate 
leaving the disc M, follows precisely the value given in Equation (I17|l calculated from the full accretor solution and not less of 
it. This is because when a fraction / of the inflowing mass accretes and the rest 1 — / remains in the disc the surface density 
E and subsequently the accretion rate increases. Within a viscous time-scale in the disc Mmax increases so that /Mmax = Md 
and not something less. Self-similar approximate solution given in Equation (|30|l describes this asymptotic stage. 
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Figure 3. Mass accretion rate in units of critical mass flow rate Mc at which partial accretion starts. The solid line shows the analytical 
expression of mass accretion rate given in Equation l|17| l. Initially M > Mc and so all mass accretes (/ = !). The dashed lines showing 
JV/niax (see Eqn. (|36j) and dotted lines showing Mt in units of Ale coincide at this stage and they both catch up with the analytical 
expression soon after the initial conditions are forgotten. When the mass inflow rate drops below Mc only a fraction / = 0.01 is allowed 
to accrete onto the central object using the boundary condition given in Equation II38I I. Although M« initially drops abruptly as a result 
of this boundary condition, it catches up with the analytical value given in Equation l|17| l in a viscous time-scale. This is because the 
density at the inner disc increases in response as 99% of the inflowing mass is retained in the disk i.e. Afmax tends to increase to M^/f 
so that its multiplication with / still gives M^. 



In order to see the behavior of the system when the boundary condition changes we should change the boundary condition 
more rapidly than the viscous time-scale in the disc. With this motivation we have changed / from 1 to 0.01 abruptly at a 
critical mass flow rate Mc at which centrifugal barrier sets in. In re al systems, of course, we expect that / will change more 
smoothly with the fastness parameter (see e.g. Eksi fc Kutlu 20ld ). The result of this simulation is shown in Figure [3] We 
observe that with the start of the partial accretion regime the accretion rate onto the star deviates from its trend given by 
Equation (|17|l . This is only for a brief epoch because the rest of the matter retained in the disc increases Mmax to M^/ f so 
that M» = /Afmax asymptotes to the value given by Equation (I17|l . The transient stage following the immediate aftermath of 
the transition to partial accretion is not modeled by the approximate self-similar solution given in Equation (|30|) as self-similar 
solutions describe the asymptotic stage when the system reaches quasi-equilibrium. 



5 DISCUSSION 

By using a recently suggested prescription we have constructed an approximate self-similar solution, given in Equation (I25p . 
describing the evolution of a partially accreting disc. This motivated another approximate solution given in Equation H30|) 
which has a constant finite inner radius. We have numerically shown that this solution describes the asymptotic stage achieved 
after a few viscous time-scales just like the exact solutions but describes the continuous range of inner boundary conditions 
between full accretion and full propeller. The approximate solution describes the inner region of the disk even better than the 

exact solution in the full accretion stage (/ = 1). 

The partial accretion might lead to a rapid decline stage ( Ibragimov fc Poutanen 200^) during an outburst of a neutron 



star soft X-ray transient when the accretion rate drops below a critical value at which the inner radius moves beyond the 
corotation radius an d the centrifugal barrier allows only a fraction of the inflowing matter from the off midplane of the disc 
l Ekgi fc Kutlulboiol ). Our results imply that partial accretion stage can change the accretion luminosity only if the decay 



timescale is short so that the inner boundary condition changes abruptly. Otherwise the mass accumulation at the inner 
boundary increases the accretion rate onto the central object to the same value it would have if there was no partial accretion 
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within the viscous time-scale. As the viscous time scale in the disc increases with the size of the disc it leads to the expectation 
that partial accretion can lead to a rapid decline stage in systems allowing small size discs i.e. transient systems with short 
orbital periods like accreting millisecond pulsars. 

The self-si milar solutio ns considered here are more appropriate for fr ee discs e.g. supernova debris discs around young 
neutron stars ([Michel fc Dessler ,1981,. ; .Alpar 2OO1I: IChatteriee et aLlboOOl ) rather than the tidally trunca ted discs in binary 
systems. A debris disc was discovered by Wang et al.l 120061 ) arouiid an anomalous X-ray pulsar (AXP; see Mereghettil (l2008h 
for a review) . Centra l compact objects (see de LucalboOSl . for a review) like Cas A might also harbor fallback discs ((Alpai||200l|). 
Tagieva et all (|2003l ) mentioned that fallback disc models o f AXPs shoul d incor porate accretion and propeller mechanisms 



working simultaneously. Partial accretion was proposed by lEksi fc Alpai 1I2OO3I ) to explain the present day luminosities of 
AXPs if they have descended from 10 ms spin periods by interacting with fallback discs. Our results imply for these systems 
that even if some of these systems are in a partial accreti on regime, the ir long term luminosity trend is determined by the 
mass flow rate derived from the full accretion solution of Pringlel (1974). Mass flow rate variations due to disc instabilities 
would result in modiflcations of this accretion rate if these changes occur sufBciently rapid. 
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